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Abstract 

We demonstrate that the integration of the recently developed dynamic mode decomposition 
(DMD) with a multi-resolution analysis allows for a decomposition method capable of robustly 
separating complex systems into a hierarchy of multi-resolution time-scale components. A 
one-level separation allows for background (low-rank) and foreground (sparse) separation of 
dynamical data, or robust principal component analysis. The multi-resolution dynamic mode 
decomposition is capable of characterizing nonlinear dynamical systems in an equation-free 
manner by recursively decomposing the state of the system into low-rank terms whose tem¬ 
poral coefficients in time are known. DMD modes with temporal frequencies near the origin 
(zero-modes) are interpreted as background (low-rank) portions of the given dynamics, and 
the terms with temporal frequencies bounded away from the origin are their sparse counter¬ 
parts. The multi-resolution dynamic mode decomposition (mrDMD) method is demonstrated 
on several examples involving multi-scale dynamical data, showing excellent decomposition re¬ 
sults, including sifting the El Nino mode from ocean temperature data. It is further applied 
to decompose a video data set into separate objects moving at different rates against a slowly 
varying background. These examples show that the decomposition is an effective dynamical 
systems tool for data-driven discovery. 


1 Introduction 

Modeling of multi-scale systems, both in time and space, pervade modern theoretical and com¬ 
putational efforts across the engineering, biological and physical sciences. Driving innovations are 
methods and algorithms that circumvent the significant challenges in efficiently connecting micro¬ 
scale to macro-scale effects that are separated potentially by orders of magnitude spatially and/or 
temporally. Wavelet-based methods and/or windowed Fourier Transforms are ideally structured to 
perform such multi-resolution analyses (MRA) as they systematically remove temporal or spatial 
features by a process of recursive refinement of sampling from the data of interest mw- Typically, 
MRA is performed in either space or time, but not both simultaneously. We propose integrating the 
concept of MRA with the recently developed Dynamic Mode Decomposition (DMD) [H [5l |6l [7J1319] , 
a technique that produces low-dimensional spatio-temporal modes. The proposed multi-resolution 
DMD (mrDMD) is shown to naturally integrate space and time so that the multi-scale spatio- 
temporal features are readily separated and approximate dynamical models constructed. 
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The origins of the DMD method, which arose from pioneering work connecting the Koopman 
operator to dynamical systems theory mm, are associated with the fluid dynamics commu¬ 
nity and the modeling of complex flows mm- Its growing success stems from the fact that it is 
an equation-free, data-driven method [3] capable of providing accurate assessments of the spatio- 
temporal coherent structures in a complex system, or short-time future estimates, thus potentially 
allowing for control protocols to be enacted simply from data sampling. The mathematical archi¬ 
tecture advocated here is an alternative to the equation-free, multi-scale modeling method proposed 
by Kevrekidis and co-workers [Ellla]. More broadly, DMD has quickly gained popularity since it 
provides information about nonlinear dynamical systems. DMD analysis can be considered to be a 
numerical approximation to Koopman spectral analysis mm, and it is in this sense that DMD is 
applicable to nonlinear systems. In fact, the terms DMD mode and Koopman mode are often used 
interchangeably in the fluids literature. 

At its core, the DMD method can be thought of as an ideal combination of the Proper Orthog¬ 
onal Decomposition (POD), a spatial dimensionality-reduction technique, with Fourier Transforms 
in time. More precisely, the DMD method produces a least-square regression to a linear dynamical 
system over the range of data collection. In this work, our goal is to integrate the DMD decom¬ 
position with key concepts from wavelet theory and MRA. Specifically, the DMD method takes 
snapshots of an underlying dynamical system to construct its decomposition. However, the fre¬ 
quency and duration (sampling window) of the data collection process can be adapted, much as in 
wavelet theory, to sift out information at different scales. Indeed, an iterative refinement of pro¬ 
gressively shorter snapshot sampling windows and recursive extraction of DMD modes from slow to 
increasingly fast time scales allows for the mrDMD. Moreover, it also allows for improved analytic 
predictions of the short-time future state of the system which is of critical importance for feedback 
control, for instance. Critical innovations demonstrated here are the ability of the mrDMD to han¬ 
dle transient phenomenon and moving (translating/rotating) structures in data, both weaknesses 
of SVD-based techniques. 

The paper is outlined as follows: In Sec. the basic DMD theory is outlined with an emphasis 
on its low-rank approximation of data. This is followed in Sec.j^by the development of the mrDMD 
structure and algorithm used in the subsequent applications, described in Sec. The applications 
exhibit the application to dynamical systems as well as more broadly to other applications such as 
video analysis. The paper is concluded in Sec. with an overview and outlook of the method. 


2 Dynamic Mode Decomposition 

The DMD method provides a spatio-temporal decomposition of data into a set of dynamic modes 
that are derived from snapshots or measurements of a given system in time. The mathematics 
underlying the extraction of dynamic information from time-resolved snapshots is closely related 
to the idea of the Arnoldi algorithm [3, one of the workhorses of fast computational solvers. The 
data collection process involves two (integer) parameters: 

N = number of spatial measurements per time snapshot, (la) 

M = number of snapshots taken in time. (lb) 

Originally the algorithm was designed to collect data at regularly spaced intervals of time. However, 
new innovations, and a more general definition of the DMD, allow for both sparse spatial [15j and 
temporal m collection of data as well as irregularly spaced collection times To illustrate the 
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algorithm, we consider regularly spaced sampling in time. The data collection times are given by: 


( 2 ) 


^m+l 


— tm + At 


where the collection time starts at ti and ends at and the interval between data collection times 
is At. In the mrDMD method, the total number of snapshots will vary as the algorithm extracts 
multi-timescale spatio-temporal structures. This will be the central focus of the next section. 

The data snapshots are arranged into an iV x M matrix 

(3) X = [x(ti) x(t 2 ) x(t 3 ) • • • x(tM)] 


where the vector x are the N measurements of the state variable of the system of interest at the data 
collection points. The objective is to mine the data matrix X for important dynamical information. 
For the purposes of the DMD method, the following matrix is also defined: 

(4) Xj^ = [x(tj) ■ ■ ■ x(tfc)]. 

Thus this matrix includes columns j through k of the original data matrix. 

The DMD method approximates the modes of the so-called Koopman operator. The Koop- 
man operator is a linear, infinite-dimensional operator that represents nonlinear, possibly infinite¬ 
dimensional, dynamics without linearization mm, and it is the adjoint of the Perron-Frobenius 
operator. The method can be viewed as computing, from the experimental data, the eigenvalues and 
eigenvectors (low-dimensional modes) of a linear model that approximates the underlying dynam¬ 
ics, even if the dynamics are nonlinear. Since the model is assumed to be linear, the decomposition 
gives the growth rates and frequencies associated with each mode. If the underlying model is linear, 
then the DMD method recovers the leading eigenvalues and eigenvectors normally computed using 
standard solution methods for linear differential equations. 

The DMD involves approximating the eigendecomposition of the best-fit linear operator A that 
relates a state Xj at time tj to the state Xj+i at the next timestep: 


(5) Xj+i«Axj. 

If Eq. ([^ holds exactly and the data is generated by a linear system, then each column of Xf^“^ 
is an element of a Krylov subspace: 


( 6 ) 


rM-1 


= [xi Axi A^xi 


, M-2 




Eq. (§ may be written in matrix form as 

(7) 


Xf 


AXf-\ 


where the operator A is chosen to minimize the Frobenius norm of ||X^ — AX^'^^^Hi?. In other 
words, the operator A advances each snapshot column in a single timestep. At, resulting in 

the future snapshot columns in X^. 

In practice, when the state dimension N is large, the matrix A may be intractable to analyze 
directly. Instead, DMD circumvents the eigendecomposition of A by considering a rank-reduced 
representation in terms of a POD-projected matrix A. The DMD algorithm proceeds as follows [9]: 
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1. First, take the SVD of 'XJf ^[17]: 

( 8 ) 


xf= usv 


where * denotes the conjugate transpose, U G ^ g (^kxk y g Here 

K is the rank of the reduced SVD approximation to The left singular vectors U are 

POD modes. 

The SVD reduction in Q could also be exploited at this stage in the algorithm to perform 
a low-rank truncation of the data. Specifically, if low-dimensional structure is present in 
the data, the singular values of S will decrease sharply to zero with perhaps only a limited 
number of dominant modes. A principled way to truncate noisy data would be to use the 
recent hard-thresholding algorithm of Gavish and Donoho [T5] . 

2. Next, compute A, the K x K projection of the full matrix A onto POD modes: 

A = XfVS-^U* 

(9) A = U*AU = U*Xf 

3. Compute the eigendecomposition of A: 

(10) AW = WA, 

where columns of W are eigenvectors and A is a diagonal matrix containing the corresponding 
eigenvalues Afe. 

4. Finally, we may reconstruct eigendecomposition of A from W and A. In particular, the 
eigenvalues of A are given by A and the eigenvectors of A (DMD modes) are given by 
columns of ’®': 

(11) = xf vs-^w. 


Note that Eq. (11) from differs from the formula ^ = UW from [5, although these will tend 
to converge if and X^ have the same column spaces. 

With the low-rank approximations of both the eigenvalues and eigenvectors in hand, the pro¬ 
jected future solution can be constructed for all time in the future. By first rewriting for convenience 
LVk = ln(Afe)/At , then the approximate solution at all future times, XDMD(i), is given by 


K 


( 12 ) 


3 bk{0)ilJk{$,) exp{u}kt) = ’®'diag(exp(a;t))b 


k=l 


where ^ are the spatial coordinates, &fc(0) is the initial amplitude of each mode, is the ma¬ 
trix whose columns are the eigenvectors ipk, diag(a;<) is a diagonal matrix whose entries are the 
eigenvalues exp(a;fct), and b is a vector of the coefficients b^. 

An alternative interpretation of ( |l^ is that it represents the least-square fit, or regression, of 
a linear dynamical system dx^Mo/di = Ax^md to the data sampled. In particular, the matrix 
A constructed is such that ||x(t) — XDMD(t)|l is minimized. In the context of the multi-resolution 
analysis that follows, each level of the multi-scale decomposition produces a linear dynamical system, 
or matrix A, for the time-scale under consideration. 
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It only remains to compute the initial coefficient values bk(0). If we consider the initial snapshot 
(xi) at time ti = 0, let’s say, then (121 gives xi = ’®'b. This generically is not a square matrix so 
that its solution 


(13) 


b = ^+xi 


can be found using a pseudo-inverse. Indeed, ’4'+ denotes the Moore-Penrose pseudo-inverse that 
can be accessed in MATLAB via the pinv command. The pseudo-inverse is equivalent to finding 
the best solution b the in the least-squares (best fit) sense. This is equivalent to how DMD modes 
were derived originally. 

Overall then, the DMD algorithm presented here takes advantage of low dimensionality in the 
data in order to make a low-rank approximation of the linear mapping that best approximates the 
nonlinear dynamics of the data collected for the system. Once this is done, a prediction of the future 
state of the system is achieved for all time. Unlike the POD-Galerkin method, which requires solving 
a low-rank set of dynamical quantities to predict the future state, no additional work is required 


for the future state prediction outside of plugging in the desired future time into (12). Thus the 


advantages of DMD revolve around the fact that (i) it is an equation-free architecture, and (ii) a 
future state prediction is known for any time t (of course, provided the DMD approximation holds). 

More broadly, the DMD method was shown to be a highly successful method for foreground/background 
subtraction in video feeds [TS]. Indeed, the DMD method is a novel, dynamical systems base method 
for performing a robust Principal Components Analysis (RPCA) of data streams |19j . Importantly, 
the DMD-based RPCA performs the low-rank/sparse matrix separation 3-4 orders of magnitude 
faster than standard ii optimization methods m- RPCA is extremely important for handling data 
outliers and/or corrupt data matrices. In the multi-resolution version of DMD, the RPCA can be 
used to effectively remove outliers at each level of decomposition. 


3 Multi-Resolution Dynamic Mode Decomposition 

The mrDMD is inspired by the observation that the slow- and fast-modes can be separated for such 
applications as foreground/background subtraction in video feeds [H]. The mrDMD recursively 
removes low-frequency, or slowly-varying, content from a given collection of snapshots. Typically, 
the number of snapshots M are chosen so that the DMD modes provide an approximately full rank 
approximation of the dynamics observed. Thus M is chosen so that all high- and low-frequency 
content is present. In the mrDMD, M is originally chosen in the same way so that an approximate 
full rank approximation can be accomplished. However, from this initial pass through the data, the 
slowest mi modes are removed, and the domain is divided into two segments with M/2 snapshots 
each. DMD is once again performed on each M/2 snapshot sequences. Again the slowest m 2 modes 
are removed and the algorithm is continued until a desired termination. 


Mathematically, the mrDMD separates the DMD approximate solution (12) in the first pass as 
follows: 


M 


M 


(14) ^ bk{0)ipi^\$.) exp{uJkt) 

k—1 k—1 k—mi-\-l 

(slow modes) (fast modes) 


where the represent the DMD modes computed from the full M snapshots. 
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The first sum in this expression (14) represents the slow-mode dynamics whereas the second 


sum is everything else. Thus the second sum can be computed to yield the matrix: 


M 


(15) 


Xm/2 — E bk{0)ip^k\^) exp(wfct). 




The DMD analysis outlined in the previous section can now be performed once again on the data 
matrix 'X.m/ 2 - However, the matrix 'X.m /2 is now separated into two matrices 

(16) = + 

where the first matrix contains the first M/2 snapshots and the second matrix contains the re- 
maining M/2 snapshots. The m 2 slow-DMD modes at this level are given hy ipl. \ where they are 
computed separately in the first of second interval of snapshots. 

The iteration process works by recursively removing slow frequency components and building 
the new matrices 'Km/ 2 : '^m/ 4 : ' ‘ ‘ until a desired/prescribed multi-resolution decomposition 

has been achieved. The approximate DMD solution can then be constructed as follows: 


(17) 


TTLi m2 

x„,DMD(t) = ^6i,^’V’i^\C)exp(a;^^^t)-l-^6fVr(^)exp(a;f^t) 

k^l k^l 

ni3 

exp(uji^h) + ■■■ 

k^l 


where at the evaluation time t, the correct modes from the sampling window are selected at each level 

(k) (k) 

of the decomposition. Specifically, the t// and oj/ are the DMD modes and DMD eigenvalues at 

(k) 

the fcth level of decomposition, the b/ ^ are the initial projections of the data onto the time interval 
of interest, and the are the number of slow-modes retained at each level. The advantage of 
this method is readily apparent: different spatio-temporal DMD modes are used to represent key 
multi-resolution features. Thus there is not a single set of modes that dominates the SVD and 
potentially marginalizes features at other time scales. 

Figure [T] illustrates the multi-resolution DMD process pictorially. In the figure, a three-level 
decomposition is performed with the slowest scale represented in blue (eigenvalues and snapshots), 
the mid-scale in red and the fast scale in green. The connection to multi-resolution wavelet analysis 
is also evident from the bottom panels as one can see that the mrDMD method successively pulls 
out time-frequency information in a principled way. 

As a final remark, the sampling strategy and algorithm discussed here (See Fig.[^ can be easily 
modified since only the slow modes at each decomposition level need to be accurately computed. 
Thus one can modify the algorithm so as to sample a fixed number, for instance M, data snapshots 
in each sampling window. The value of M need not be large as only the slow modes need to be 
resolved. Thus the sampling rate would increase as the decomposition proceeds from one level to 
the next. This assures that the lowest levels of the mrDMD are not highly-sampled since the cost 
of the SVD would be greatly increased by such a fine sampling rate. 


3.1 Formal mrDMD Expansion 

The solution 0 can be made more precise. Specifically, one must account for the number of levels 
(L) of the decomposition, the number of time bins (J) for each level, and the number of modes 
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Figure 1: Representation of the multi-resolution dynamic mode decomposition where successive 
sampling of the data, initially with M snapshots and decreasing by a factor of two at each resolution 
level, is shown (top figures). The DMD spectrum is shown in the middle panel where there are 
mi (blue dots) slow-dynamic modes at the slowest level, m 2 (red) modes at the next level and m 3 
(green) modes at the fastest time-scale shown. The shaded region represents the modes that are 
removed at that level. The bottom panels shows the wavelet-like time-frequency decomposition of 
the data color coded with the snapshots and DMD spectral representations. 


retained at each level (mL). This can be easily seen in Fig. Thus the solution is parametrized 
by the following three indices: 

t' = 1, 2, • • • , L number of decomposition levels (18a) 

j = 1, 2, • • • , J number time bins per level (J = 2^^“^^) (18b) 

A: = 1,2, • • • , mz, number of modes extracted at level L. (18c) 

To formally define the series solution for x„„dmd(A), we propose the following indicator function 
(19) hAt) = [l elsewlere'^ with j = and J = 

which is only non-zero in the interval, or time bin, associated with the value of j. The parameter 
i denotes the level of the decomposition. 
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Figure 2: Illustration of the mrDMD mode decomposition and hierarchy. Represented are the 
modes tpkH^) ^^.d their position in the decomposition structure. The triplet of integer values, £,j 
and k, uniquely express the time level, bin and mode of the decomposition. 


The three indices and indicator function (19) give the mrDMD solution expansion 


( 20 ) 


>(^) = E E E exp(4^’^'^t). 


{=1j=ik=i 


This is a concise definition of the mrDMD solution that includes the information on the level, time 
bin location and number of modes extracted. Figure demonstrates the mrDMD decomposition 
in terms of the solution (20). In particular, each mode is represented in its respective time bin and 


level. An alternative interpretation of this solution is that it yields the least-square fit, at each level 
£ of the decomposition, to the linear dynamical system 


( 21 ) 


dx^^d) 

dt 


— _4(t,i)jj.(f,i) 


where the matrix captures the dynamics in a given time bin j at level £. 

The indicator function acts as sifting function for each time bin. Interestingly, this 

function acts as the Gabor window of a windowed Fourier transform [3]. Since our sampling bin 
has a hard cut-off of the time series, it may introduce some artificial high-frequency oscillations. 
Time-series analysis, and wavelets in particular, introduce various functional forms that can be used 
in an advantageous way. Thus thinking more broadly, one can imagine using wavelet functions for 
the sifting operation, thus allowing the time function fe,j(t) to take the form of one of the many 
potential wavelet basis, i.e. Haar, Daubechies, Mexican Hat, etc. This will be considered in future 


work. For the present, we simply use the sifting function introduced in (19) 






















Figure 3: Comparison of the true modes {ijjk) to DMD (-(/ifc) and mrDMD (ipl,) modes. The true, 
mrDMD and DMD modes are showed in columns from left to right. The inset shows the relative 
error of the mrDMD = W'ipj — and DMD {Rj = W'ipj — il’jW) in comparison to the 

true mode. Temporal dynamics are also compared, with true dynamics shown in red dots, DMD 
reconstructed dynamics shown in blue solid line and mrDMD reconstructed dynamics shown in 
black solid line. 

4 Application of Method 

The mrDMD developed in the last section is implemented here on three example data sets: the 
first being an exemplar of a video stream, the second coming from atmospheric-ocean data, and 
the third an exemplar of data with moving objects. This final example is especially important 
as it renders many SVD-based decompositions useless. For the first case, the implementation is 
compared against the standard DMD algorithm, highlighting the ability of the mrDMD to correctly 
capture multi-scale phenomena. It should be noted that in previous work with DMD, which was 
not framed in the mrDMD architecture advocated here, a level-1 decomposition was effectively used 
to separate the foreground from the background in a video m- Indeed, this work motivates the 
generalization to the multi-resolution analysis. 
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Figure 4: Comparison of snapshots from the original video of Fig.j^and reconstruction from DMD 
and mrDMD. The true dynamics are showed in the first row, while the mrDMD reconstructed 
dynamics are showed in the second row and the DMD reconstructed dynamics are showed in the 
third row. The time labels on the bottom correspond to when each snapshot was taken. The 
mrDMD reconstruction matches with the original video very well while the DMD reconstruction 
shows inconsistencies. Indeed, the mrDMD method can easily handle signals that turn on and off 
in time while the DMD method cannot. 


4.1 Spatio-temporal filtering of video 

The first example we consider for application of the mrDMD is illustrated in Fig. For this 
example, four spatio-temporal modes are combined into a single data set. Specifically, we combine 
the four modes shown in the left panels with the time dynamics given in the right panels. The four 
modes used to construct the true solution are represented by ipj{x,y) for j = 1,2,3 and 4. Their 
corresponding time dynamics are given by aj (t). Thus the true solution is expressed by 


( 22 ) 


= ■ 

i=i 


Both the DMD (represented by and the modes ipj of (121 with j = 1,2,3 and 4) and mrDMD 

(represented by x,„,dmd and the of (17) where fc = 1 and i = 1,2,3) attempt to reconstruct 

X. 

Figure shows the true modes along with the approximating mrDMD and DMD modes. The 
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Figure 5: Eigenvalues from level-1 to level-3 of the mrDMD decomposition applied to the data 
generated in Fig. This is a specific case of the mrDMD abstraction shown in Fig.[^ The shaded 
circle shows the threshold radius used for background mode subtraction. The inset illustrates the 
time-frequency position of each sampling window. Also illustrated are the mode selection and 
labeling process. 


inset in the mrDMD and DMD modes shows the £2 difference between the true modes and the 
approximating modes of the respective decomposition, i.e. Rj^i = W-ipj — ^-nd Rj = W-ipj —i’jW 

respectively. The mrDMD modes are almost identical to the true modes while the DMD modes 
show an error in mixing the modes. The error in the DMD modes is clearly illustrated in the time 
dynamics. Specifically, the DMD does not correctly capture the on-off dynamics of modes three 
and four so that ||x — Xomd|| 0(1). In contrast, the mrDMD is able to easily capture the correct 
time dynamics with its windowing decomposition. Indeed, the agreement between the mrDMD and 
exact solution is remarkable, i.e. ||x — x,„dmd|| ‘C 1. 

To further illustrate the accuracy of the mrDMD decomposition in comparison to the DMD 
decomposition, consider Fig. This shows the true solution x at the time snap shots of t = 
25,50, 75,100 and 125. The mixing of the four modes of (22) is nicely captured by the mrDMD 
while the DMD fails to capture key features at these time points. Indeed, at time t = 75 in 
particular, the agreement is quite poor between the DMD and true solution. 


As a final part of the analysis of the video sequence represented by x(t) in (22), we can explicitly 
perform the multi-resolution decomposition depicted in Fig.j^ In this case, the dominant eigenvalue 
are kept at each level of the decomposition. The bottom centered box is the level I decomposition 
showing that there exists a single mode in the mode removal region (shaded blue inset). This is 
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Figure 6 : Application of mrDMD on sea surface temperature data from 1990 to 2010. The left panel 
illustrates the process for a 4-level decomposition. At each level, the slowest modes are extracted. 
At a given level, the zero-mode component has a period T = oo. Illustrated are the (a) Level-1 
mrDMD mode with period T = oo, and (b) Level-1 mrDMD mode with period T = 52 weeks. 
These two modes are in the gray mode removal region of Fig. Further in the decomposition we 
can extract (c) Level-4 mrDMD mode of 1997 with period T = oo and (d) Level-4 mrDMD mode 
of 1999 with period T = oo. Mode (c) clearly shows the El Nino mode of interest that develops 
in the central and east-central equatorial Pacific. The El Nino mode was not present in 1999 as is 
clear from mode (d). Data source: NOAA_OI_SST_V2 data provided by the NOAA/OAR/ESRL 
PSD, Boulder, Colorado, USA, from their Web site at http://www.esrl.noaa.gov/psd/ 

mode At the next level, there are two modes (complex conjugate pairs) in the mode removal 

region (shaded pink inset), and (l^ft panel) and and (right panel). These 

modes represent the slow oscillatory dynamics driven by 02 (t). And in the final level, an additional 
pair of eigenvalues are present in the first three frames of the decomposition that correspond to the 
dynamics of modes 03 ( 1 ) (frames 1 and 2) and 04 ( 1 ) (frame 3) with modes and Thus 

the multi-resolution strategy outlined in Fig. holds. Moreover, the intuition derived from this 
wavelet-like strategy is remarkably effective in decomposing the multi-resolution, spatio-temporal 
dynamics of x. 
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4.2 Multi-scale time separation of complex system: El Nino, Southern 
Oscillation 

The example of the previous subsection was contrived to demonstrate the ability of the mrDMD 
to separate the spatio-temporal modes of Fig. In this example, we use the mrDMD on a more 
realistic data set. Specifically, we consider global surface temperature data over the ocean. The data 
is open source from the NOAA/OAR/ESRL PSD, Boulder, Colorado, USA. The NOAA_OI_SST_V2 
data set considered can be downloaded at http://www.esrl.noaa.gov/psd/ The data spans a 20 year 
period from 1990 to 2010. 

Figure shows the results of the mrDMD algorithm. Specifically, a 4-level decomposition is 
performed with the slow spatio-temporal modes pulled at each level as suggested in Fig. The 
zero mode in each window is the DC component, or period infinity T = oo mode. At the Hrst level 
of the decomposition, two modes are extracted: the zero mode (T = oo) depicted in Fig. |^a) and a 
yearly cycle {T = 52 weeks) shown in Fig. j^b). The yearly cycle is the slowest mode extracted at 
level 1. Note that the zero mode component of level 1, is just the average ocean temperature 

over the entire 20-year data set. 

We can continue the mrDMD analysis through to the fourth-level. At the fourth level, the 
approximate zero mode (period T = oo) of the sampling window extracts physically interesting 
results. In particular, the data-driven method of the mrDMD discovers the 1997 El Nifio mode 
generated from the well-known El Nifio, Souther Oscillation (ENSO). Indeed, 1997 is known to have 
been a strong El Nino year, as verified by its strong modal signature in the 4th level decomposition 
of the mrDMD. In contrast, the same sampling window shifted down to 1999 produces no El Nino 
mode, which is in keeping with known ocean patterns that year. El Nifio is the warm phase of the 
ENSO cycle and is associated with a band of warm ocean water that develops in the central and 
east-central equatorial Pacihc (between approximately the International Date Line and 120°W), 
including off the Pacific coast of South America. The mrDMD mode clearly shows this band of 
warm ocean water as a spatio-temporal mode in 1997 in Fig. |^c). ENSO refers to the cycle of 
warm and cold temperatures, as measured by sea surface temperature, SST, of the tropical central 
and eastern Pacific Ocean. El Nifio is accompanied by high air pressure in the western Pacific and 
low air pressure in the eastern Pacihc. The cool phase of ENSO is called La Nina with SST in the 
eastern Pacihc below average and air pressures high in the eastern and low in western Pacihc. The 
ENSO cycle, both El Nifio and La Nina, causes global changes of both temperatures and rainfall. 

These results could not have been produced with DMD unless the correct sampling windows 
were chosen ahead of time, thus requiring a supervised learning step not required by the mrDMD. 
Further, even in such a case, the previous slow modes, such as those of level-1 demonstrated in 
Fig.[^ a) and (b) would pollute the data at the level of investigation. Thus the mrDMD provides a 
principled, algorithmic approach that is capable of data-driven discovery in complex systems such 
as ocean/atmospheric data. 

4.3 Translating and/or Rotating Strnctnres 

The final application of the mrDMD is on an example that is notoriously difficult for SVD-based 
methods to characterize, namely, when translational and/or rotational structures are present in the 
data snapshots, i.e. continuous or discrete invariances. Indeed, such invariances completely under¬ 
mine our ability to compute low-rank embeddings of the data as driven by correlated structures, or 
POD/PCA modes. This has been the Achilles heel of many SVD based methods, thus requiring in 
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Figure 7: Application of mrDMD on moving objects separation. Two modes (top panel labeled 
“Fast Mode” and “Slow Mode”) are combined in the data snapshots (top right panel). The “Fast 
Mode” moves from left to right as indicated by the arrow with a speed of lOv while the “Slow Mode” 
moves from bottom to top at the speed v. We take v = 1/40 without loss of generality. Thus the 
fast and slow mode speeds are approximately an order of magnitude different. In the mrDMD 
decomposition (left bottom panel), the “Slow Mode” is extracted at level 2 as represented by the 
panels (a) and (b). The “Fast Mode” is extracted at level 13 as represented at three representative 
panels (c)-(e). Although there is some shadow (residual) of the slow mode on the fast mode and 
vice-versa, the mrDMD is remarkably effective in extracting the correct modes. Moreover, the level 
at which they are extracted can allow for a reconstruction of the velocity and direction. To our 
knowledge, this is the best performance achieved to date with an SVD-based method with multiple 
time-scale objects. 


applications such as PCA-based face recognition (eigenfaces) [H], well-cropped and centered faces 
for reasonable performance, i.e. translation and rotation are removed in an expensive pre-processing 
procedure. 

In a dynamical systems setting, a simple traveling wave will appear to be a high-dimensional 
object in POD space despite the fact that it is only constructed from two modes, one associated 
with translational invariance. For dynamical cases exhibiting translation and/or rotation, Rowley 
and Marsden [52] formulated one of the only mathematical strategies to date to extract the low¬ 
dimensional embeddings. In particular, they developed a template-matching technique to first 
remove the invariance before applying the SVD decomposition. Although effective, it is not suited 
for cases where a myriad of objects and timescales are present in the data. 

The mrDMD is well suited to handle invariances such as translation and rotation. Consider 
once again the case of a simple traveling wave. The standard DMD decomposition applied to the 
traveling wave would result in a solution approximation requiring many DMD modes due to the 
slow fall off of the singular values Q in step 1 of the DMD algorithm. Further, the eigenvalues 
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ujk in (121 would also be bounded away from the origin as there is no background mode for such 
translating data. In the mrDMD architecture, the fact that the eigenvalues are bounded away 
from the origin (and typically 0(1)) in the initial snapshot window Xm would simply allow the 
mrDMD method to ignore the traveling wave at the first level of the mrDMD. Once the domain is 
now divided into two for the next level of analysis, the traveling wave is now effectively moving at 
half the speed in this new domain, i.e. the eigenvalues have migrated towards the origin and the 
traveling wave is now re-evaluated. The recursive procedure would eventually produce a sampling 
window where the traveling wave looks sufficiently stationary and low-rank so as to be extracted 
in the multi-resolution analysis. The level at which it is extracted also characterizes the speed of 
the traveling wave. Specifically, the higher the level in the decomposition where the traveling wave 
is extracted, the faster its speed. 

Figure demonstrates the application of the mrDMD method on a simple example in which 
there are two moving objects, one moving at a slow velocity and another moving at a high velocity. 
Specifically, the example constructed results from the dynamical sequence 


(23) x = tpi{x - 10vt,y) +'ijj 2 {x,y - vt) 

where the two modes (Gaussians of the form = exp(—< 7 ( 2 ; — xq)^ — a(y — yo)^) with a = 0.1 
and (xo,yo) = (—18,20) and (xo,yo) = (-20,-9) for the fast and slow modes respectively) used 
to construct the true solution are represented by 'iljj{x,y) for j = 1,2. Note that the first mode is 
translating from left to right at speed lOu whereas the second mode is translating from bottom to 
top at speed v. Without loss of generality, v can be set to any value. It is chosen to be u = 1/40 
for the domain and Gaussians considered. In this case, the straightforward template matching 
procedure would fail due to the two distinct time-scales of the objects. As shown in the figure, the 
mrDMD is capable of pulling out the two modes at level 2 (slow) and 13 (fast) of the analysis. It 
is at these levels that the translating objects look close to stationary (Iw/jj ^ 1) in the mrDMD 
analysis. The level 2 modes correspond to the slow-moving object whereas the level 13 modes are 
associated with the fast object. Although there is a residual in both the extracted slow and fast 
modes, it does a reasonable job in extracting the fast and slow modes. To our knowledge, this is 
the only SVD-based method capable of such unsupervised performance. 

The motivation for this example is quite clear when considering video processing and surveil¬ 
lance. In particular, for many applications in this arena, background subtraction is only one step 
in the overall video assessment. Identification of the foreground objects is the next, and ultimately 
more important, task. Gonsider the example of a video of a street scene in which there is a pedes¬ 
trian walking (slow translation) along with a vehicle driving down the road (fast translation). Not 
only would we want to separate the foreground from background, we would also want to separate 
the pedestrian from the vehicle. More precisely, we would want to make a separate video and 
identification of the slow and fast objects in the video. The results of Fig.j^show the mrDMD to be 
well-suited and effective in this task. Moreover, one can also envision augmenting the algorithm with 
a recursive application of mrDMD to improve performance much like in the foreground/background 
separation work on DMD m- 


5 Discussion and Outlook 

Data-driven strategies for analyzing complex systems are of growing interest in the mathematical 
sciences. Indeed, methods capable of providing principled decompositions of data arising in multi¬ 
scale, spatio-temporal systems are key enabling strategies for many applications in the engineering, 
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physical and biological sciences. The multi-resolution dynamic mode decomposition advocated 
in this work capitalizes on recent innovations in equation-free modeling via the dynamic mode 
decomposition. It leverages these ideas by integrating them with concepts from wavelet theory 
and multi-resolution analysis. By construction, the method provides a principled reconstruction of 
multi-resolution, spatio-temporal data sets. The methods effectiveness is demonstrated on several 
example data sets, highlighting its ability to extract critical information and enact data-driven 
discovery protocols. 

The DMD method provides a spatio-temporal decomposition of data into a set of dynamic 
modes that are derived from snap shots or measurements of a given system in time. The DMD 
method approximates the modes of the so-called Koopman operator. The Koopman operator is a 
linear, infinite-dimensional operator that represents nonlinear, infinite-dimensional dynamics with¬ 
out linearization, and is the adjoint of the Perron-Frobenius operator. The method can be viewed as 
computing, from the experimental data, the eigenvalues and eigenvectors (low-dimensional modes) 
of a linear model that approximates the underlying dynamics, even if the dynamics is nonlinear. 
By interpreting the DMD eigenvalues as corresponding to prescribed time scale dynamics, one can 
extract spatio-temporal structures recursively for shorter and shorter sampling windows. Thus the 
slow-modes are removed first and the data is filtered for analysis of its higher frequency content. 
This recursive sampling structure is demonstrated to be effective in allowing for a reconstruction 
of several example data sets. 

One can envision a number of innovations to augment the proposed mrDMD strategy. Many 
are particularly attractive for applications across the engineering, physical and biological sciences. 
Indeed, the impact that DMD is having on complex fluid flows is already known and has already 
been mentioned in the introduction. Fields like neuroscience, which are rich in multi-scale, complex 
dynamics are also ideal candidates for exploration using the mrDMD infrastructure as DMD has 
already had recent demonstrated success in this arena [23]. The outlook of these techniques is 
highlighted here: 

Compressive sampling: First is the ability to leverage tools from compressive sampling m to 
facilitate the collection of considerably fewer measurements, resulting in the same multi-resolution 
dynamic mode decomposition, but with considerably fewer measurements. This reduction in the 
number of measurements may have a broad impact in situations where data acquisition is expensive 
and/or prohibitive. In particular, we envision these tools being used in particle image velocimetry 
(PIV) to reduce the data transfer requirements for each snapshot in time, increasing the maximum 
temporal sampling rate. Other applications include ocean and atmospheric monitoring, where 
individual sensors are expensive. Even if full-state measurements are available, the proposed method 
of compressed DMD will be computationally advantageous in situations where there is low-rank 
structure in the high-dimensional data. 

DMD Control: A second important direction revolves around recent innovations of the DMD 
with control [24] . which is capable of disambiguating between the underlying dynamics and the 
effects of actuation, resulting in accurate input-output models. The method is data-driven in that 
it does not require knowledge of the underlying governing equations, only snapshots of state and 
actuation data from historical, experimental, or black-box simulations. One can envision developing 
such input-output models at the various levels of spatio-temporal encoding in the data, i.e. multi¬ 
resolution input-output models can be assessed and constructed from a given system. Both these 
research directions, and many others, highlight the potential strength of the model in pushing 
forward principled, equation-free strategies for the analysis of complex multi-scale systems. 
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Machine Learning: Modern tools of statistical analysis and dimensionality-reduction have become 
the workhorses for the burgeoning field of machine learning (ML). ML techniques aim to capitalize 
on underlying low-dimensional patterns and clustering in data. In the dynamical applications 
considered here, one might exploit these patterns, or DMD modes, by building libraries of low- 
rank dynamical modes, much like is done with POD modes PSI \T7\ . Such DMD libraries 
for different dynamical regimes partner nicely with compressive sensing strategies. Additionally, 
Kernel based techniques, which are at the core of support vector machines, for instance, have 
already found successful application in the DMD architecture when considering more accurate, 
nonlinear dynamical reconstructions [28) . Maximum advantage should be taken of such techniques 
when integrating the DMD architecture in applications. 
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